Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

cicecore: correct initial condition metadata #818

Merged
merged 2 commits into from
Mar 13, 2023

Conversation

phil-blain
Copy link
Member

@phil-blain phil-blain commented Mar 8, 2023

PR checklist

  • Short (1 sentence) summary of your PR:
    Make the initial condition metadata more correct by writing it at the initial time and not marking its variables as averaged.
  • Suggest PR reviewers from list in the column to the right.
  • Please copy the PR test results link or provide a summary of testing completed below.
    I did not run the test suite, as this only changes history output and the our regression test compare the restarts, not the history output. I checked the following:
    • that the initial condition is now written with the correct date/time in the file name (the model initialization date/time)
    • that the time netCDF variable in the initial condition has value zero and no bounds attribute
    • that the time_bounds variable is not written in the initial condition
    • that the variable in the initial condition have time_rep = "instantaneous" and no cell_method attribute.
    • that the rest of the outputs variables, both for the initial condition and the simulation output, are unchanged. The history files can't be compared with cmp since the time of the run is written as a global :history attribute, but I used something like
 diff <(ncdump  history/iceh_01h.2005-01-01-03600.nc) <(ncdump  history-ic-fixed/iceh_01h.2005-01-01-03600.nc)

to verify that the only difference is this global history attribute.

  • How much do the PR code changes differ from the unmodified code?
    • bit for bit (expect for changes described above)
    • different at roundoff level
    • more substantial
  • Does this PR create or have dependencies on Icepack or any other models?
    • Yes
    • No
  • Does this PR update the Icepack submodule? If so, the Icepack submodule must point to a hash on Icepack's main branch.
    • Yes
    • No
  • Does this PR add any new test cases?
    • Yes
    • No
  • Is the documentation being updated? ("Documentation" includes information on the wiki or in the .rst files from doc/source/, which are used to create the online technical docs at https://readthedocs.org/projects/cice-consortium-cice/. A test build of the technical docs will be performed as part of the PR testing.)
    • Yes
    • No, does the documentation need to be updated at a later time?
      • Yes
      • No
  • Please provide any additional information or relevant details below:

I reported these two issues initially in #567 (comment):

[...] this leads to the initial condition being written after the first time step...

Another minor thing I noticed: when using hist_avg = .true., the initial condition is also written as being "averaged" from "time = 0" to "time = after the first time step", which is kind of weird...

I can remove the changes to the standalone driver if we do not agree on that, but in my opinion it does not break anything and makes the history output more correct so I felt it was better to just do it there also.

Full commit messages follow.


[PATCH 1/2] ice_history_write: fix initial condition metadata under 'hist_avg'

When writing averaged history outputs (hist_avg=.true.), this setting
also affects the initial condition. Even if the actual data variables
written to the initial condition are not averaged (they are taken
more or less directly from the restart or the hard-coded defaults,
modulo aggregation over categories), their attributes ('cell_method' and
'time_rep') imply they are averaged, and the 'bound' attribute of the
'time' variable refers to the 'time_bounds' variable.

Make the metadata of the initial condition more correct by:
- not writing the 'time_bounds' (and the corresponding 'd2' dimension)
- not writing the 'bounds' attribute of the 'time' variable
- not writing the 'cell_method' attributes of each variable
- writing the 'time_rep' attribute of each variable as 'instantaneous'
instead of 'averaged'.

Do this by checking 'write_ic' at all places where we check for the
value of 'hist_avg' to write the above variables and attributes in each
of the 3 IO backends (binary, netcdf, pio2).
[PATCH 2/2] drivers/{nemo_concepts,standalone}: write initial condition at initial time

In CICE_InitMod::cice_init, we call ice_calendar::advance_timestep
before writing the initial condition, such that the 'time' variable in
the initial condition is not zero; it has a value of 1*dt (the model
time step). The initial condition filename also reflects this, since
'msec' (model seconds) also has a value of 1*dt and is used in
ice_history_shared::construct_filename. This leads to the initial
condition filename not corresponding to the model initialization
date/time but rather 1*dt later.

Fix that by calling 'accum_hist' to write the initial condition _before_
calling 'advance_timestep'.

Copy link
Contributor

@dabail10 dabail10 left a comment

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Looks good.

@apcraig
Copy link
Contributor

apcraig commented Mar 9, 2023

Just to clarify, is the history ic file data identical compared to baseline? Did just the filename and some attributes change? Or is the file written at a different part of the initialization sequence or even prior to advancing now?

@phil-blain
Copy link
Member Author

Thanks for your questions, it made me take a closer look at the standalone code.

The data is now written prior to advancing the calendar by calling advance_timestesp, so that is why time changes. It is still written as part of the initialization sequence, prior to stepping the model, but yes, at a different place in cice_init.

The diff below shows the whole context in between the addition and the deletion for the standalone driver:

       ! Initialize shortwave components using swdn from previous timestep
       ! if restarting. These components will be scaled to current forcing
       ! in prep_radiation.
       if (trim(runtype) == 'continue' .or. restart) &
          call init_shortwave    ! initialize radiative transfer

+      if (write_ic) call accum_hist(dt) ! write initial conditions
+
 ! tcraig, use advance_timestep here
 !      istep  = istep  + 1    ! update time step counters
 !      istep1 = istep1 + 1
 !      time = time + dt       ! determine the time and date
 !      call calendar(time)    ! at the end of the first timestep
       call advance_timestep()

    !--------------------------------------------------------------------
    ! coupler communication or forcing data initialization
    !--------------------------------------------------------------------

       call init_forcing_atmo    ! initialize atmospheric forcing (standalone)

       if (tr_fsd .and. wave_spec) call get_wave_spec ! wave spectrum in ice
       call get_forcing_atmo     ! atmospheric forcing from data
       call get_forcing_ocn(dt)  ! ocean forcing from data

       ! snow aging lookup table initialization
       if (tr_snow) then         ! advanced snow physics
          call icepack_init_snow()
          call icepack_warnings_flush(nu_diag)
          if (icepack_warnings_aborted()) call abort_ice(error_message=subname, &
             file=__FILE__, line=__LINE__)
          if (snw_aging_table(1:4) /= 'test') then
             call init_snowtable()
          endif
       endif

       ! isotopes
       if (tr_iso)     call fiso_default                 ! default values

       ! aerosols
       ! if (tr_aero)  call faero_data                   ! data file
       ! if (tr_zaero) call fzaero_data                  ! data file (gx1)
       if (tr_aero .or. tr_zaero)  call faero_default    ! default values
       if (skl_bgc .or. z_tracers) call get_forcing_bgc  ! biogeochemistry
       if (z_tracers) call get_atm_bgc                   ! biogeochemistry

       if (runtype == 'initial' .and. .not. restart) &
          call init_shortwave    ! initialize radiative transfer using current swdn

       call init_flux_atm        ! initialize atmosphere fluxes sent to coupler
       call init_flux_ocn        ! initialize ocean fluxes sent to coupler

-      if (write_ic) call accum_hist(dt) ! write initial conditions
-
       if (my_task == master_task) then
          call ice_memusage_print(nu_diag,subname//':end')
       endif

       end subroutine cice_init

So, I think that since the forcing is read after the initial condition is written (with my changes), if any output variable that correspond to forcing fields (e.g. f_uatm, etc) are written to the history, then yes the data in the initial condition will change with this PR: now it will probably be zero (since the forcing has not been read). This is not super nice, I guess.

It's not the case in our NEMO driver since all forcing is passed to CICE by NEMO.

So an alternative, instead of moving the accum_hist call up, would be to move the call to advance_timestep down, after writing the initial condition. I could do that.

Apart from that, the changes are mostly in attributes, with this change, in the initial condition :

  • the time_bounds variable (and the corresponding d2 dimension) is gone
  • the time variable has value zero, before, it had value dt
  • the time variable does not have the 'bounds' attribute
  • the variables do not have a cell_method attribute
  • the time_rep attribute of each variable is 'instantaneous' instead of 'averaged'.
  • has the model initialization date/time as filename (instead of that+dt)

@phil-blain
Copy link
Member Author

So, I think that since the forcing is read after the initial condition is written (with my changes), if any output variable that correspond to forcing fields (e.g. f_uatm, etc) are written to the history, then yes the data in the initial condition will change with this PR: now it will probably be zero (since the forcing has not been read).

I've just checked that's indeed the case.

@phil-blain
Copy link
Member Author

So, with this diff instead, the filename and attributes are OK (time=0), etc., and the data stays the same as it is now:

diff --git a/./cicecore/drivers/standalone/cice/CICE_InitMod.F90 b/./cicecore/drivers/standalone/cice/CICE_InitMod.F90
index 8de05a1..c94b431 100644
--- a/./cicecore/drivers/standalone/cice/CICE_InitMod.F90
+++ b/./cicecore/drivers/standalone/cice/CICE_InitMod.F90
@@ -192,27 +192,20 @@ subroutine cice_init

       if (tr_aero .or. tr_zaero) call faero_optics !initialize aerosol optical
                                                    !property tables

       ! Initialize shortwave components using swdn from previous timestep
       ! if restarting. These components will be scaled to current forcing
       ! in prep_radiation.
       if (trim(runtype) == 'continue' .or. restart) &
          call init_shortwave    ! initialize radiative transfer

-! tcraig, use advance_timestep here
-!      istep  = istep  + 1    ! update time step counters
-!      istep1 = istep1 + 1
-!      time = time + dt       ! determine the time and date
-!      call calendar(time)    ! at the end of the first timestep
-      call advance_timestep()
-
    !--------------------------------------------------------------------
    ! coupler communication or forcing data initialization
    !--------------------------------------------------------------------

       call init_forcing_atmo    ! initialize atmospheric forcing (standalone)

       if (tr_fsd .and. wave_spec) call get_wave_spec ! wave spectrum in ice
       call get_forcing_atmo     ! atmospheric forcing from data
       call get_forcing_ocn(dt)  ! ocean forcing from data

@@ -238,20 +231,28 @@ subroutine cice_init
       if (z_tracers) call get_atm_bgc                   ! biogeochemistry

       if (runtype == 'initial' .and. .not. restart) &
          call init_shortwave    ! initialize radiative transfer using current swdn

       call init_flux_atm        ! initialize atmosphere fluxes sent to coupler
       call init_flux_ocn        ! initialize ocean fluxes sent to coupler

       if (write_ic) call accum_hist(dt) ! write initial conditions

+! tcraig, use advance_timestep here
+!      istep  = istep  + 1    ! update time step counters
+!      istep1 = istep1 + 1
+!      time = time + dt       ! determine the time and date
+!      call calendar(time)    ! at the end of the first timestep
+      call advance_timestep()
+
+
       if (my_task == master_task) then
          call ice_memusage_print(nu_diag,subname//':end')
       endif

       end subroutine cice_init

I think this is a better approach, Let me know what you think and I can update the PR.

@apcraig
Copy link
Contributor

apcraig commented Mar 9, 2023

Thanks @phil-blain. I'm a little surprised that shifting the call of advance timestep later in the sequence doesn't change the model results. The model time will be different when some of those init methods are called. It could be they are just setting up the run, not actually initializing the data based on model time. In that case, that's great. But we should probably run a full test suite just to make sure the answers don't change. I can help do that if you like, just let me know.

@phil-blain
Copy link
Member Author

Yes, I was planning on running a full suite. I'll fire that up right now.

@eclare108213
Copy link
Contributor

Thanks for fixing this, @phil-blain! When writing of the initial condition was added many years ago, the whole initialization process was pretty complicated and different in various coupled models, and I was mainly interested in seeing the initial ice state, not so much the forcing. So this is an improvement. I definitely think it should be tested not only in standalone mode but also the various coupled configurations, to make sure they're all BFB in their solutions, if not the initial condition file output itself.

@phil-blain
Copy link
Member Author

OK, so as hinted yesterday during our call, moving the call to advance_timestep to the end of cice_init is not b4b with master, and worse, it leads to several restart tests in the base_suite to fail. So this is not as easy as I thought :P

So I see two options:

  1. I drop the changes to the standalone driver from this PR. Then it stays written with (initial date+dt) as the filename and time=dt, as it is currently.
  2. I keep the PR as is (just moving the call to accum_hist before advance_timestep). This corrects the filename and time variable, but then forcing variables in the initial condition will be zero (or whatever their default value is). Thinking about it more, I think this is better than the current situation. With the current situation, the data in the initial condition is not all at the same date, since the ice state (from the restart) is at time=0 but the forcing is at time=dt, so there is a bit of an incoherence.

Let me know what you think.

@daveh150
Copy link
Contributor

daveh150 commented Mar 10, 2023 via email

@phil-blain
Copy link
Member Author

There is also option 3:

  1. initialize the forcing before calling advance_timestep, write the initial condition, advance the timestep and re-initialize the forcing with the new date. This would allow the forcing to be written in the initial condition at the initial date, but should not change model answers. But this is a little bit invasive so I think I still prefer option 2 as its simpler.

@apcraig
Copy link
Contributor

apcraig commented Mar 10, 2023

I think we don't want to change model answers unless we do a careful refactor of the standalone sequencing with respect to forcing, time advancement, and model timestepping. I'm open to doing that, but it's a bigger effort.

I think option 2 is a good compromise at this point, but am open. I think @eclare108213 should weigh in too.

@eclare108213
Copy link
Contributor

Is option 2 what is currently coded? This looks okay to me.

@phil-blain
Copy link
Member Author

Thanks @eclare108213. Yes, option 2 is what is currently in this PR. I'll run a full suite just to be extra careful, I don't expect answers to change.

When writing averaged history outputs (hist_avg=.true.), this setting
also affects the initial condition. Even if the actual data variables
written to the initial condition are not averaged (they are taken
more or less directly from the restart or the hard-coded defaults,
modulo aggregation over categories), their attributes ('cell_method' and
'time_rep') imply they are averaged, and the 'bound' attribute of the
'time' variable refers to the 'time_bounds' variable.

Make the metadata of the initial condition more correct by:
- not writing the 'time_bounds' variable (and the corresponding 'd2' dimension)
- not writing the 'bounds' attribute of the 'time' variable
- not writing the 'cell_method' attributes of each variable
- writing the 'time_rep' attribute of each variable as 'instantaneous'
instead of 'averaged'.

Do this by checking 'write_ic' at all places where we check for the
value of 'hist_avg' to write the above variables and attributes in each
of the 3 IO backends (binary, netcdf, pio2).
…l time

In CICE_InitMod::cice_init, we call ice_calendar::advance_timestep
before writing the initial condition, such that the 'time' variable in
the initial condition is not zero; it has a value of 1*dt (the model
time step). The initial condition filename also reflects this, since
'msec' (model seconds) also has a value of 1*dt and is used in
ice_history_shared::construct_filename. This leads to the initial
condition filename not corresponding to the model initialization
date/time but rather 1*dt later.

Since we call 'accum_hist' after initializing the forcing, any forcing
field written to the initial condition has values corresponding to
msec=dt, whereas the ice state corresponds to msec=0, leading to an
inconsistency.

Fix that by calling 'accum_hist' to write the initial condition _before_
calling 'advance_timestep'. Since we now call 'accum_hist' before
initializing the forcing, any forcing field written to the initial
condition have its default, hard-coded value, instead of its value at
time=dt. An improvement would be to read the forcing at time=dt, write
the initial condition, advance the time step, and read the forcing
again, but let's not complicate things too much for now.
@phil-blain
Copy link
Member Author

Base suite is b4b:

349 measured results of 349 total results
349 of 349 tests PASSED
0 of 349 tests PENDING
0 of 349 tests MISSING data
0 of 349 tests FAILED

I've just pushed an updated version, only commit messages are changed, range-diff:

1:  610c8a9 ! 1:  9f4504d ice_history_write: fix initial condition metadata under 'hist_avg'
    @@ Commit message
         'time' variable refers to the 'time_bounds' variable.

         Make the metadata of the initial condition more correct by:
    -    - not writing the 'time_bounds' (and the corresponding 'd2' dimension)
    +    - not writing the 'time_bounds' variable (and the corresponding 'd2' dimension)
         - not writing the 'bounds' attribute of the 'time' variable
         - not writing the 'cell_method' attributes of each variable
         - writing the 'time_rep' attribute of each variable as 'instantaneous'
2:  25b7c2b ! 2:  d72cdf9 drivers/{nemo_concepts,standalone}: write initial condition at initial time
    @@ Commit message
         condition filename not corresponding to the model initialization
         date/time but rather 1*dt later.

    +    Since we call 'accum_hist' after initializing the forcing, any forcing
    +    field written to the initial condition has values corresponding to
    +    msec=dt, whereas the ice state corresponds to msec=0, leading to an
    +    inconsistency.
    +
         Fix that by calling 'accum_hist' to write the initial condition _before_
    -    calling 'advance_timestep'.
    +    calling 'advance_timestep'. Since we now call 'accum_hist' before
    +    initializing the forcing, any forcing field written to the initial
    +    condition have its default, hard-coded value, instead of its value at
    +    time=dt. An improvement would be to read the forcing at time=dt, write
    +    the initial condition, advance the time step, and read the forcing
    +    again, but let's not complicate things too much for now.

      ## cicecore/drivers/direct/nemo_concepts/CICE_InitMod.F90 ##
     @@ cicecore/drivers/direct/nemo_concepts/CICE_InitMod.F90: subroutine cice_in

@apcraig apcraig merged commit 9a55ad9 into CICE-Consortium:main Mar 13, 2023
phil-blain added a commit to phil-blain/CICE that referenced this pull request Mar 16, 2023
* ice_history_write: fix initial condition metadata under 'hist_avg'

When writing averaged history outputs (hist_avg=.true.), this setting
also affects the initial condition. Even if the actual data variables
written to the initial condition are not averaged (they are taken
more or less directly from the restart or the hard-coded defaults,
modulo aggregation over categories), their attributes ('cell_method' and
'time_rep') imply they are averaged, and the 'bound' attribute of the
'time' variable refers to the 'time_bounds' variable.

Make the metadata of the initial condition more correct by:
- not writing the 'time_bounds' variable (and the corresponding 'd2' dimension)
- not writing the 'bounds' attribute of the 'time' variable
- not writing the 'cell_method' attributes of each variable
- writing the 'time_rep' attribute of each variable as 'instantaneous'
instead of 'averaged'.

Do this by checking 'write_ic' at all places where we check for the
value of 'hist_avg' to write the above variables and attributes in each
of the 3 IO backends (binary, netcdf, pio2).

* drivers/{nemo_concepts,standalone}: write initial condition at initial time

In CICE_InitMod::cice_init, we call ice_calendar::advance_timestep
before writing the initial condition, such that the 'time' variable in
the initial condition is not zero; it has a value of 1*dt (the model
time step). The initial condition filename also reflects this, since
'msec' (model seconds) also has a value of 1*dt and is used in
ice_history_shared::construct_filename. This leads to the initial
condition filename not corresponding to the model initialization
date/time but rather 1*dt later.

Since we call 'accum_hist' after initializing the forcing, any forcing
field written to the initial condition has values corresponding to
msec=dt, whereas the ice state corresponds to msec=0, leading to an
inconsistency.

Fix that by calling 'accum_hist' to write the initial condition _before_
calling 'advance_timestep'. Since we now call 'accum_hist' before
initializing the forcing, any forcing field written to the initial
condition have its default, hard-coded value, instead of its value at
time=dt. An improvement would be to read the forcing at time=dt, write
the initial condition, advance the time step, and read the forcing
again, but let's not complicate things too much for now.

(cherry picked from commit 9a55ad9)

Cherry-pick notes:

There were some conflicts in io_{netcdf,pio2}/ice_history_write.F90 due
to 26d917a (Fix QC test, fix bug in history time axis, fix history
output averaging for timestep output (CICE-Consortium#624), 2021-08-19), since that
commit enabled averaging over multiple time steps and thus removed the
assumption that histfreq(ns) = '1' means hist_avg = .false.. I also had
to add additional changes to 'ice_history_write'  in the conditions that
are checked before writing the NetCDF attributes since we do not yet
have the 'ice_write_hist_attrs' subroutine added in 078aab4 (Merge
cgridDEV branch including C grid implementation and other fixes (CICE-Consortium#715),
2022-05-10).

Closes: https://gitlab.science.gc.ca/concepts/CICE/issues/19
phil-blain added a commit to phil-blain/CICE that referenced this pull request Mar 16, 2023
* ice_history_write: fix initial condition metadata under 'hist_avg'

When writing averaged history outputs (hist_avg=.true.), this setting
also affects the initial condition. Even if the actual data variables
written to the initial condition are not averaged (they are taken
more or less directly from the restart or the hard-coded defaults,
modulo aggregation over categories), their attributes ('cell_method' and
'time_rep') imply they are averaged, and the 'bound' attribute of the
'time' variable refers to the 'time_bounds' variable.

Make the metadata of the initial condition more correct by:
- not writing the 'time_bounds' variable (and the corresponding 'd2' dimension)
- not writing the 'bounds' attribute of the 'time' variable
- not writing the 'cell_method' attributes of each variable
- writing the 'time_rep' attribute of each variable as 'instantaneous'
instead of 'averaged'.

Do this by checking 'write_ic' at all places where we check for the
value of 'hist_avg' to write the above variables and attributes in each
of the 3 IO backends (binary, netcdf, pio2).

* drivers/{nemo_concepts,standalone}: write initial condition at initial time

In CICE_InitMod::cice_init, we call ice_calendar::advance_timestep
before writing the initial condition, such that the 'time' variable in
the initial condition is not zero; it has a value of 1*dt (the model
time step). The initial condition filename also reflects this, since
'msec' (model seconds) also has a value of 1*dt and is used in
ice_history_shared::construct_filename. This leads to the initial
condition filename not corresponding to the model initialization
date/time but rather 1*dt later.

Since we call 'accum_hist' after initializing the forcing, any forcing
field written to the initial condition has values corresponding to
msec=dt, whereas the ice state corresponds to msec=0, leading to an
inconsistency.

Fix that by calling 'accum_hist' to write the initial condition _before_
calling 'advance_timestep'. Since we now call 'accum_hist' before
initializing the forcing, any forcing field written to the initial
condition have its default, hard-coded value, instead of its value at
time=dt. An improvement would be to read the forcing at time=dt, write
the initial condition, advance the time step, and read the forcing
again, but let's not complicate things too much for now.

(cherry picked from commit 9a55ad9)

Cherry-pick notes:

There were some conflicts in io_{netcdf,pio2}/ice_history_write.F90 due
to 26d917a (Fix QC test, fix bug in history time axis, fix history
output averaging for timestep output (CICE-Consortium#624), 2021-08-19), since that
commit enabled averaging over multiple time steps and thus removed the
assumption that histfreq(ns) = '1' means hist_avg = .false.. I also had
to add additional changes to 'ice_history_write'  in the conditions that
are checked before writing the NetCDF attributes since we do not yet
have the 'ice_write_hist_attrs' subroutine added in 078aab4 (Merge
cgridDEV branch including C grid implementation and other fixes (CICE-Consortium#715),
2022-05-10).
phil-blain added a commit to phil-blain/CICE that referenced this pull request Mar 16, 2023
phil-blain added a commit to phil-blain/CICE that referenced this pull request Feb 2, 2024
* ice_history_write: fix initial condition metadata under 'hist_avg'

When writing averaged history outputs (hist_avg=.true.), this setting
also affects the initial condition. Even if the actual data variables
written to the initial condition are not averaged (they are taken
more or less directly from the restart or the hard-coded defaults,
modulo aggregation over categories), their attributes ('cell_method' and
'time_rep') imply they are averaged, and the 'bound' attribute of the
'time' variable refers to the 'time_bounds' variable.

Make the metadata of the initial condition more correct by:
- not writing the 'time_bounds' variable (and the corresponding 'd2' dimension)
- not writing the 'bounds' attribute of the 'time' variable
- not writing the 'cell_method' attributes of each variable
- writing the 'time_rep' attribute of each variable as 'instantaneous'
instead of 'averaged'.

Do this by checking 'write_ic' at all places where we check for the
value of 'hist_avg' to write the above variables and attributes in each
of the 3 IO backends (binary, netcdf, pio2).

* drivers/{nemo_concepts,standalone}: write initial condition at initial time

In CICE_InitMod::cice_init, we call ice_calendar::advance_timestep
before writing the initial condition, such that the 'time' variable in
the initial condition is not zero; it has a value of 1*dt (the model
time step). The initial condition filename also reflects this, since
'msec' (model seconds) also has a value of 1*dt and is used in
ice_history_shared::construct_filename. This leads to the initial
condition filename not corresponding to the model initialization
date/time but rather 1*dt later.

Since we call 'accum_hist' after initializing the forcing, any forcing
field written to the initial condition has values corresponding to
msec=dt, whereas the ice state corresponds to msec=0, leading to an
inconsistency.

Fix that by calling 'accum_hist' to write the initial condition _before_
calling 'advance_timestep'. Since we now call 'accum_hist' before
initializing the forcing, any forcing field written to the initial
condition have its default, hard-coded value, instead of its value at
time=dt. An improvement would be to read the forcing at time=dt, write
the initial condition, advance the time step, and read the forcing
again, but let's not complicate things too much for now.

(cherry picked from commit 9a55ad9)

Cherry-pick notes:

There were some conflicts in io_{netcdf,pio2}/ice_history_write.F90 due
to 26d917a (Fix QC test, fix bug in history time axis, fix history
output averaging for timestep output (CICE-Consortium#624), 2021-08-19), since that
commit enabled averaging over multiple time steps and thus removed the
assumption that histfreq(ns) = '1' means hist_avg = .false.. I also had
to add additional changes to 'ice_history_write'  in the conditions that
are checked before writing the NetCDF attributes since we do not yet
have the 'ice_write_hist_attrs' subroutine added in 078aab4 (Merge
cgridDEV branch including C grid implementation and other fixes (CICE-Consortium#715),
2022-05-10).
phil-blain added a commit to phil-blain/CICE that referenced this pull request Feb 2, 2024
phil-blain added a commit to phil-blain/CICE that referenced this pull request Feb 9, 2024
* ice_history_write: fix initial condition metadata under 'hist_avg'

When writing averaged history outputs (hist_avg=.true.), this setting
also affects the initial condition. Even if the actual data variables
written to the initial condition are not averaged (they are taken
more or less directly from the restart or the hard-coded defaults,
modulo aggregation over categories), their attributes ('cell_method' and
'time_rep') imply they are averaged, and the 'bound' attribute of the
'time' variable refers to the 'time_bounds' variable.

Make the metadata of the initial condition more correct by:
- not writing the 'time_bounds' variable (and the corresponding 'd2' dimension)
- not writing the 'bounds' attribute of the 'time' variable
- not writing the 'cell_method' attributes of each variable
- writing the 'time_rep' attribute of each variable as 'instantaneous'
instead of 'averaged'.

Do this by checking 'write_ic' at all places where we check for the
value of 'hist_avg' to write the above variables and attributes in each
of the 3 IO backends (binary, netcdf, pio2).

* drivers/{nemo_concepts,standalone}: write initial condition at initial time

In CICE_InitMod::cice_init, we call ice_calendar::advance_timestep
before writing the initial condition, such that the 'time' variable in
the initial condition is not zero; it has a value of 1*dt (the model
time step). The initial condition filename also reflects this, since
'msec' (model seconds) also has a value of 1*dt and is used in
ice_history_shared::construct_filename. This leads to the initial
condition filename not corresponding to the model initialization
date/time but rather 1*dt later.

Since we call 'accum_hist' after initializing the forcing, any forcing
field written to the initial condition has values corresponding to
msec=dt, whereas the ice state corresponds to msec=0, leading to an
inconsistency.

Fix that by calling 'accum_hist' to write the initial condition _before_
calling 'advance_timestep'. Since we now call 'accum_hist' before
initializing the forcing, any forcing field written to the initial
condition have its default, hard-coded value, instead of its value at
time=dt. An improvement would be to read the forcing at time=dt, write
the initial condition, advance the time step, and read the forcing
again, but let's not complicate things too much for now.

(cherry picked from commit 9a55ad9)

Cherry-pick notes:

There were some conflicts in io_{netcdf,pio2}/ice_history_write.F90 due
to 26d917a (Fix QC test, fix bug in history time axis, fix history
output averaging for timestep output (CICE-Consortium#624), 2021-08-19), since that
commit enabled averaging over multiple time steps and thus removed the
assumption that histfreq(ns) = '1' means hist_avg = .false.. I also had
to add additional changes to 'ice_history_write'  in the conditions that
are checked before writing the NetCDF attributes since we do not yet
have the 'ice_write_hist_attrs' subroutine added in 078aab4 (Merge
cgridDEV branch including C grid implementation and other fixes (CICE-Consortium#715),
2022-05-10).

Rebase onto CICE6.4.1 notes:

There were some conflicts in io_{netcdf,pio2}/ice_history_write.F90 due
to the same commits mentioned above. They are both included in
CICE6.4.1, but the original cherry-picked commit (9a55ad9 (cicecore:
correct initial condition metadata (CICE-Consortium#818), 2023-03-13)) is not.
The correct resolution was thus to:
- keep ours version for calls to ice_write_hist_attrs
- cherry-pick from 9a55ad9 the changes to ice_write_hist_attrs
- in io_pio2, change 'if (hist_avg)' to 'if (hist_avg .and. .not.
  write_ic)'
phil-blain added a commit to phil-blain/CICE that referenced this pull request Feb 9, 2024
@phil-blain phil-blain deleted the write-ic-at-init-date branch February 19, 2024 20:16
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Projects
None yet
Development

Successfully merging this pull request may close these issues.

5 participants